Computation and visualization of photonic quasicrystal spectra via Block's theorem 
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Previous methods for determining photonic quasicrystal (PQC) spectra have relied on the use 
of large supercells to compute the eigenfrequencies and/or local density of states (LDOS). In this 
manuscript, we present a method by which the energy spectrum and the eigenstates of a PQC can be 
obtained by solving Maxwell's equations in higher dimensions for any PQC defined by the standard 
cut-and-project construction, to which a generalization of Bloch's theorem applies. In addition, 
we demonstrate how one can compute band structures with defect states in the higher-dimensional 
superspace with no additional computational cost. As a proof of concept, these general ideas are 
demonstrated for the simple case of one-dimensional quasicrystals, which can also be solved by 
simple transfer-matrix techniques. 

PACS numbers: 



I. INTRODUCTION 

We propose a computational method to solve for 
the spectra and eigenstates of quasicrystalline elec- 
tromagnetic structures by directly solving a periodic 
eigenproblem in a higher-dimensional lattice. Such 
photonic quasicrystals (PQCs) have a number of 
unique properties compa red to ordinary p eriodic struc- 
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in two or three dimensions where they can have greater 
rotational symmetry and therefore offer some hope of 
complete photonic band gaps with lower index con- 
traslP^ ^^ I ^^ I than the roughly 2:1 contrast currently re- 
quired for periodic structures. However, the study of 
two- and three-dimensional photonic quasicrystals has 
been hampered by the computational difficulty of mod- 
eling aperiodic structures, which has previously required 
large "supercell" calculations that capture only a portion 
of the infinite aperiodic lattice. Our method, in contrast, 
captures the entire infinite aperiodic structure in a sin- 
gle higher-dimensional unit cell, and we believe that this 
approach will ultimately be much more computationally 
tractable for two- and three-dimensional quasicrystals. 
The idea that many quasicrystals can be constructed by 
an irr ational s lice of a higher-dimensional lattice is well 
known ,'^^1^^'^ and in fact is the most common f ormula- 
tion of quasicrystals in two and three dimensions, 
but the possibility of direct numerical calculations within 
the higher-dimensional space seems to have been little 
explored outside of some tight-binding calculations in 
quantum systems.^^ As a proof of concept, we demon- 
strate a first implementation of the technique applied to 
one-dimensional quasicrystals, such as the well known Fi- 
bonacci structure. Not only can we reproduce the spec- 
trum from transfer-matrix calculations, but we also show 
that the higher-dimensional picture provides an interest- 
ing way to visualize the eigenmodes and compute defect 
states in the infinite aperiodic structure. 

There have been several previous numerical approaches 



to simulating quasicrystal structures in electromagnetism 
and quantum mechanics. In one dimension, a typical 
quasicrystal is an aperiodic sequence of two or more 
materials, determined either by a slice of a higher- 
dimensional lattic^^ or by some "string concatenation" 
rule.^^ In either case, efficient 2x2 transfer-matrix meth- 
ods are available that allow one to quickly compute the 
transmission spectra and density of stat es for supercells 
consisting of many thousands of layers. '^^'^ Two- and 
three-dimensional quasicrystals are almost always de- 
fined as an irrational slice (i.e., incommensurate Miller 
indices) of a higher-dimensional lattice; for example, the 
famous Penrose tiling can be viewed as a two-dimensional 
slice of a five-dimensional cubic lattice or of a four- 
dimensional root lattice A/^?^ In such cases, supercell 
computations of a finite portion of the infinite aperi- 
odic structure (or a rational approximant thereof) ^^'^ 
require slower numerical methods, most commonly finite- 
difference time-doma in (F DTD) simulations^^ or 
planewave expansions .'^^'^ Unfortunately, these methods 
become very expensive for large supercells, nearly pro- 
hibitively so for three-dimensional quasic rystals — there 
have been experiments for 3D PQCs,'^^^ but as yet few 
theoretical predictions.^^ With FDTD methods, for ex- 
ample, the PQC local density of states is typically inte- 
grated in Mon te- Car lo fashion via random sources or ini- 
tial conditions'^'^'^^'but many simulations are required to 
sample all possible modes in a large supercell. Also, the 
finite domain of a supercell becomes even more significant 
in higher dimensions where a tractable superc ell is nec- 
essarily smaller, as there can be localized stateJ^m^^^'^ 
whose presence is dependent on the particular region of 
the PQC considered. Our method of computing the spec- 
trum directly in the higher-dimensional unit cell, on the 
other hand, requires no supercell to capture the infinite 
aperiodic structure — it uniformly samples (up to a finite 
resolution) every possible supercell of the infinite qua- 
sicrystal, rather than any particular subsection. The 
influence of finite-resolution on the convergence of the 



2 



spectrum can be systematically understood: one is not 
"missing" any part of the quasicrystal, so much as re- 
solving the entire quasicrystal with lower resolution. 

The structure of this paper is as follows: in Sec. [n| 
we review the "cut-and-project" method for defining a 
PQC as a slice of a higher-dimensional lattice, followed 
in Sec. |III| by a description of our computational method 
in the higher-dimensional lattice. There, we describe 
the extension of Maxwell's equations to higher dimen- 
sions and also describe its solution in terms of a higher- 
dimensional Bloch planewave expansion. As a proof of 
concept, we present a sequence of one-dimensional ex- 



amples in Sec. [IVj First, we compare results for a one- 
dimensional "Fibonacci sequence" with standard one- 
dimensional transfer-matrix techniques. Second, as men- 
tioned above, cut-and-project allows for a straightfor- 
ward way of studying defects in the the quasicrystal with 
the same computational effort as the perfect PQC, and 
this is demonstrated in the one-dimensional "Fibonacci" 
example. Finally, we demonstrate the ease with which 
one can construct and explore different quasicrystals by 
continuously varying the cut angle. 




FIG. 1: Unit cell of the Fibonacci superspace dielectric. The 
physical dielectric is obtained by taking a slice at an angle 
tan(/) = T. Black/ white are the dielectric constants of the 
structure factor material and air, chosen to be e = 4.84 and 
e — 2.56, respectively. 



II. QUASICRYSTALS VIA CUT-AND-PROJECT 

Given a periodic lattice, any lower dimensional cross- 
section of that lattice may be either periodic or quasi- 
periodic, depending upon the angle of the cross-section. 
For example, the periodic 2D cross-sections of a 3D crys- 
tal are the lattice planes, defined in crystallography by 
integer Miller indices. If the Miller indices have irrational 
ratios, on the other hand, the cross-section is aperiodic 
but still has long-range order because of the underlying 
higher-dimensional periodicity. This is what is known as 
a "cut-and-project" method of defining a quasicrystalline 
structure: as a slice of a periodic structure in a higher- 
dimensional "superspace" j ^^ l ^'^ l (For a thorough discus- 
sion of quasicrystals via cut-and-project, see Ref. [26l ) 
Cut-and-project defines a specific class of quasicrystals; 
equivalently, and more abstractly, cut-and-project corre- 
sponds to structures whose Fourier transform has sup- 
port spanned by a finite number of reciprocal basis vec- 
tors (the projection of th e reciprocal lattice vectors from 
higher dimensions) j ^^ l ^^ l This class includes most com- 
monly considered quasicrystals in two or three dimen- 
sions, including the Penrose tiling, as well as many 
one-dimensional quasicrystals including a version of the 
Fibonacci structure. 

For example, consider the Fibonacci PQC in one di- 
mension formed from two materials sa = 4.84 and £b = 
2.56 in layers of thickness A and 5, respectively, simi- 
lar to a recent experimental structure.^ The Fibonacci 
structure S is then defined by the limit n ^ oo of the 
string-concatenation rule Sn = Sn-2Sn-i with starting 
strings So = B and = A,^ generating a sequence 
BABAABABAABA • • • . In the case where B/A is the 
golden ratio r = (1 + a/5)/2, exactly the same structure 



can be generated by a slice of a two-dimensional lattice 
as depicted in Fig. [l]^^ The slice is at an angle (j) with 
an irrational slope tan^ = 1/r, and the unit cell of the 
2D lattice is an A x A square at an angle ^ in a square 
lattice with period {A + B) sin cj) = a. Because the slope 
is irrational, the offset /intercept of the slice is unimpor- 
tant: any slice at an angle (j) intercepts the unit cell at 
infinitely many points, filling it densely. 

For thickness ratios B/A ^ r, the Fibonacci structure 
cannot be constructed by cut-and-project, and in general 
string-concatenation rules can produce a different range 
of structures (such as the Thue-Morse PQC)^than cut- 
and-project. This is partly a question of definition — 
some authors reserve th e term "quasicrystal" for cut- 
and-project structures.^ In any case, cut-and-project in- 
cludes a wide variety of aperiodic structures, including 
most of the structures that have been proposed in two 
or three dimensions (where they can be designed to have 
n-fold rotational symmetry for any n), and are the class 
of quasicrystals that we consider in this paper. 

In general, let d < 3 be the number of physical di- 
mensions of a quasicrystal structure generated by a d- 
dimensional "slice" of an n-dimensional periodic struc- 
ture {n > d). Denote this slice by X (the physical space) 
with coordinates x G M^, and denote the remaining n — d 
coordinates by y G W^~^ in the "unphysical" space Y (so 
that the total n-dimensional superspace is Z = X F). 
The primitive lattice vectors 'R.i G Z define the orienta- 
tion of the lattice with respect to the slice (rather than 
vice versa), with corresponding primitive reciprocal vec- 
tors Gi defined by the usual • G^- = 27r6ij^ (The 
concept of an "irrational slice" is commonly used in the 
quasicrystal literature. However, a general definition of 
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what is meant by an "irrational slice" seems difficult to 
find, and less evident in dimensions d > 2. For a more 
precise definition of "irrational slice" in general dimen- 
sions and a proof that it is dense in the unit cell, see 
Sec.|Vl) 

The physical dielectric function £:(x) is then con- 
structed by starting with a periodic dielectric function 
£:(x, y) in the superspace and evaluating it at a fixed y 
(forming the slice) . Because an irrational slice is dense in 
the unit cell of the superspace, it doesn't matter what 
value of y one chooses, as discussed below (as long as 
e is piecewise continuous). In principle, one could de- 
fine the unit cell of e in the superspace to be any arbi- 
trary n-dimensional function, but in practice it is com- 
mon to "decorate" the higher-dimension unit cell with 
extrusions of familiar d-dimensional objects.^^ More 
precisely, "cut-and-project" commonly refers to construc- 
tions where a set of lattice points within a finite window 
of the cut plane are projected onto the cut plane, and 
this is equivalent to a simple cut where objects at the 
lattice points are extruded in the y direction by the win- 
dow widthP^In particular, the extrusion window is com- 
monly an inverted projection (shadow) of the unit cell 
onto the y directions, although this is not the case for 
the Fibonacci construction of Fig. [l] 

(Note that the higher-dimensional lattice need not be 
hypercubic. For example, the Penrose tiling can be 
expressed as a two-dimensional slice of either a five- 
dimensional hypercubic lattice or of a non-orthogonal 
four-dimensional root lattice ^4.^^ For computational 
purposes, the lower the dimensionality the better.) 

III. COMPUTATIONS IN HIGHER 
DIMENSIONS 

Although the cut-and-project technique is a standard 
way to define the quasicrystal structure, previous com- 
putational studies of photonic quasicrystals then pro- 
ceeded to simulate the resulting structure only in the 
projected (d-dimensional) physical space. Instead, it is 
possible to extend Maxwell's equations into the periodic 
n-dimensional superspace, where Bloch's theorem applies 
to simplify the computation. By looking at only the 
unit cell in n dimensions one can capture the infinite d- 
dimensional quasicrystal. Our development of this tech- 
nique was inspired by earlier research on analogous elec- 
tronic quasicrystals that applied a tight-binding method 
in two dimensions to compute th e spe ctrum of a one- 
dimensional electronic quasicrystal .^^^^ 

Let us start with Maxwell's equations in the physi- 
cal space X for the quasicrystal 6:(x, y) at some fixed y 
(that is, y is viewed as a parameter, not a coordinate). 
Maxwell's equations can be written as an eigenproblem 
for the harmonic modes H(x, y)e~*'^*,l^ where again y 
appears as a parameter. 

Vx X ^^Vx X H = {co/cfH, (1) 



where Vx x denotes the curl with respect to the x coor- 
dinates. Assuming that the structure is quasicrystalline, 
i.e. that X is an irrational slice of the periodic super- 
space Z, then CO should not depend upon y.^^ The reason 
is that y only determines the offset of the "initial" slice 
of the unit cell (for x = 0), but as we reviewed above 
the slice (considered in all copies of the unit cell) fills 
the unit cell densely. Therefore, any change of y can be 
undone, to arbitrary accuracy, merely by offsetting x to 
a different copy of the unit cell. An offset of x doesn't 
change the eigenvalues cj, although of course it offsets the 
eigenfunctions H. 

The fact that uj is independent of y allows us to re- 
interpret Eq. ([1]), without actually changing anything: 
we can think of y as a coordinate rather than a pa- 
rameter, and the operator on the left-hand side as an 
operator in d-dimensional space. Note that H is still a 
three-component vector field, and Vx x is still the or- 
dinary curl operator along the x directions, so this is 
not so much a higher-dimensional version of Maxwell's 
equations as an extension of the unmodified ordinary 
Maxwell's equations into a higher-dimensional parame- 
ter space. The y coordinate appears in the operator only 
through e. Because uj is independent of y, i.e. it is 
just a number rather than a function of the coordinates, 
the equation ([T]) in higher dimensions is still an eigen- 
problem, and its spectrum of eigenvalues co is the same 
as the spectrum of the d-dimensional quasicrystal, since 
the equations are identical. The physical solution is ob- 
tained by evaluating these higher-dimensional solutions 
at a fixed y, say y = (where a different y merely cor- 
responds to an offset in x as described above). 

For a real, positive £, both the physical operator and 
the extended operator in in Eq. are Hermitian and 
positive semi-definite, leading to many important prop- 
erties such as real frequencies cj.^^ 



A. Bloch's theorem and numerics for quasicrystals 

Because the superspace eigenproblem is periodic, 
Bloch's theorem applies: the eigenfunctions H(x, y) can 
be written in the Bloch form h(z)e*^"^, where h is a pe- 
riodic function defined by its values in the unit cell, and 
k is the n-dimensional Bloch wavevector.^^ 

Here, k determines the phase relationship between H 
in different unit cells of the superspace, but it does not 
have a simple interpretation once the solution is pro- 
jected into physical space. The reason is that h, viewed 
as a function of x, is again only quasiperiodic: transla- 
tion in X "wraps" the slice into a different portion of the 
unit cell, so both h and e*^'^ change simultaneously and 
the latter phase cannot be easily distinguished. This pre- 
vents one from defining a useful phase or group velocity 
of the PQC modes. 

The key point is that Bloch's theorem reduces the 
eigenproblem to a finite domain (the n-dimensional unit 
cell), rather than the infinite domain required to describe 



4 



the quasicrystal solutions in physical space. This means 
that standard numerical methods to find the eigenvalues 
of differential operators are immediately applicable. For 
example, since the solution h is periodic, one can apply 
a planewave expansion method^^ for h: 



h(z) = ^hGe* 



iG z 



(2) 



where the summation is over all n-dimensional recipro- 
cal lattice vectors G. Because the curl operations only 
refer to the x coordinates, Vx x h is replaced by a sum- 
mation over Gx x he, where Gx denotes G projected 
into X. The resulting eigenproblem for the Fourier co- 
efficients h (once they are truncated to some wavevector 
cutoff) can be computed either by direct dense-matrix 
methods"^^ or, more efficiently, by iterative methods ex- 
ploiting fast Fourier transforms. In the present paper, 
we do the former, which is easy to implement as a proof 
of concept, but for higher-dimensional computations an 
iterative method will become necessary. 

We should also remind the reader that there is a con- 
straint Vx • H = on the eigenfunctions, in order to ex- 
clude unphysical solutions with static magnetic charges. 
In a planewave method, this leads to a trivial constraint 
(kx + Gx) • h = 0, again with k and G projected into X. 



B. The spectrum of the quasicrystal 

With a familiar eigenproblem arising from Bloch's the- 
orem, such as that of a periodic physical structure, the 
eigenvalues form a band structure: discrete bands c^n(k) 
that are continuous functions of k, with a finite number 
of bands in any given frequency range.^^ For a finite- 
resolution calculation, one obtains a finite number of 
these bands cOn with some accuracy that increases with 
resolution, but even at low resolutions the basic struc- 
ture of the low- frequency bands is readily apparent. The 
eigenvalues of the higher-dimensional quasicrystal oper- 
ator of Eq. ([T]), on the other hand, are quite different. 

The underlying mathematical reason for the discrete 
band structure of a physical periodic structure is that 
the Bloch eigen-operator for a periodic physical lattice, 
(V + zk) X ^ ( V + ik) X , is the inverse of a compact in- 
tegral operator corresponding to the Green's function, 
and hence the spectral theorem applies.!^ Among other 
things, this implies that the eigenvalues at any given k 
for a finite unit cell form a discrete increasing sequence, 
with a finite number of eigenvalues below any finite co. 
The same nice property does not hold for the operator 
extended to n dimensions, because along the y directions 
we have no derivatives, only a variation of the scalar func- 
tion e. Intuitively, this means that the fields can oscil- 
late very fast along the y directions without necessarily 
increasing co^ allowing one to have infinitely many eigen- 
functions in a finite bandwidth. More mathematically, 
an identity operator is not compact and does not satisf 
the spectral theorem,^ and since the operator of Eq. (Il 



is locally the identity along the y directions the same 
conclusion applies. This means that, when the y direc- 
tion is included as a coordinate, it is possible to get an 
infinite number of bands in a finite bandwidth at a fixed 
k. 

In fact, as we shall see below, this is precisely what 
happens, and moreover it is what must happen in order 
to reproduce the well-known properties of quasicrystal 
spectra. It has been shown that quasicrystal spectra can 
exhibit a fractal structure, with infinitely many gaps 
(of decreasing size) in a finite bandwidth, and such a 
structure could not arise from an ordinary band diagram 
with a finite number of bands in a given bandwidth. Of 
course, once the unit cell is discretized for numerical com- 
putation, the number of degrees of freedom and hence the 
number of eigenvalues is finite. However, as the resolu- 
tion is increased, not only do the maximum frequency 
and the accuracy increase as for an ordinary computa- 
tion, but also the number of bands in a given bandwidth 
increases. Thus, as the resolution is increased, more and 
more of the fractal structure of the spectrum is revealed. 



IV. ONE-DIMENSIONAL RESULTS 

As a proof of concept implementation of cut-and- 
project, we construct a Fibonacci quasicrystal in 
Sec. |IV A| using the projection method described above, 
compute the band structure as a function of the projected 
wave- vector and compare to a transfer-matrix calcula- 
tion of the same quasicrystal structure. We also demon- 
strate the field visualization enabled by the projection 
method, both in the superspace (n dimension s) as w ell 
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as in the physical space {d dimensions). In Sec. 
demonstrate how this method can accommodate systems 
with defects without additional computational costs. Fi- 
nally, we explore several one-dimensional quasicrystal 
configurations in Sec. |IV C] by varying the cut angle (j). 



A. Fibonacci quasicrystal 

We solved Eq. ([T]) numerically using a planewave ex- 
pansion in the unit cell of the 2D superspace, as described 
above, for the ID Fibonacci quasicrystal structure de- 
picted in Fig. [l] The resulting band diagram is shown 
in Fig. [2|left), along with a side-by-side comparison of 
the local density of states in Fig. [2|right) calculated us- 
ing a transfer-matrix approach with a supercell of 10^ 
layers The two calculations show excellent agreement 
in the location of the gaps, except for one or two easily- 
identified spurious bands inside some of the gaps, which 
are discussed in further detail below. The most impor- 
tant feature of Fig. [2|left) is the large number of bands 
even in the finite bandwidth uo G [0,0.4], with the num- 
ber of bands increasing proportional to the spatial res- 
olution (planewave cutoff). This is precisely the feature 



predicted abstractly above, in Sec. IIIB at a low resolu- 
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(27r/a) p(co) 

FIG. 2: Left: Frequency spectrum oo of the Fibonacci quasicrystal vs. "wave- vector" kx- The blue Unes indicate spurious 
states which arise due to finite resolution effects (see text). Right: Corresponding density of states p{(jj) computed using a 
transfer-matrix technique with a supercell of 10^ layers. Inset: Power distribution ^ \Hz\'^ (red/white = positive/zero) of a 
spurious state, the field profile of which oscillates at the Nyquist frequency. 
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FIG. 3: Enlarged view of the Fibonacci spectrum showing a 
gap with a spurious band crossing it. Insets show the mag- 
netic field \Hz\ for the spurious band at various kx — the lo- 
calization of this mode around the X-parallel edges of the 
dielectric indicate that this is a discretization artifact. 



tion, one sees only the largest gaps, and at higher resolu- 
tions further details of the fractal spectrum are revealed 
as more and more bands appear within a given band- 
width, very different from calculations for periodic phys- 
ical media. The important physical quantity is not so 
much the band structure, since k has no simple physical 
meaning as discussed previously, but rather the density 
of states formed by projecting the band structure onto 



the CO axis. In this density of states, the small number of 
spurious bands within the gaps, which arise from the dis- 
cretization as discussed below, plays no significant role: 
the density of states is dominated by the huge number 
of flat bands (going to infinity as the resolution is in- 
creased), and the addition of one or two spurious bands 
is negligible. 

The "spurious" bands that appear within some of the 
band gaps of the superspace calculation arise from the 
discretization of the dielectric interfaces parallel to the 
slice direction. Because the slice is at an irrational an- 
gle, it will never align precisely with a uniform grid, re- 
sulting in inevitable staircasing effects at the boundary. 
With ordinary electromagnetic simulations, these stair- 
casing effects can degrade the accuracy, but here the 
lack of derivatives perpendicular to the slice allows spuri- 
ous modes to appear along these staircased edges (there 
is no frequency penalty to being localized perpendicular 
to the slice). Indeed, if one looks at the field patterns for 
the spurious modes, one of which is shown in the inset of 
Fig. [2|right), one sees that the field intensity is peaked 
along the slice-parallel dielectric interfaces. Because they 
are localized to these interfaces and therefore dominated 
by the unphysical staircasing, the spurious modes behave 
quite differently from the "real" solutions and are easily 
distinguished qualitatively and quantitatively. Most im- 
portantly, as the resolution is increased, the number of 
spurious modes in a given gap does not increase like all of 
the other bands, because the thickness of the staircased 
interface region decreases proportional to the resolution. 
This makes the gaps in the band structure obvious: here, 
they are the only frequency ranges for which the number 
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FIG. 4: Plot of the magnetic field amplitude \Hz\ for a 
band-edge state taken along a slice of the two-dimensional 
superspace (in the (j) direction). Inset: Two-dimensional 
superspace field profile (red/ white/blue indicates posi- 
tive / zero / negative amplitude) . 

of eigenvalues does not increase with resolution. Equiv- 
alently, as noted above, the contribution of the spurious 
bands to the density of states is asymptotically negligible 
as resolution is increased. 

Computing the eigenmodes in the higher-dimensional 
superspace immediately suggests a revealing visualiza- 
tion technique: instead of plotting the quasiperiodic 
fields as a function of the physical coordinates x by tak- 
ing a slice, plot them in the two-dimensional superspace. 
This has the advantage of revealing the entire infinite 
aperiodic field pattern in a single finite plot.l^ One such 
plot was already used above, to aid in understanding the 
spurious modes localized at staircased interfaces. A typ- 
ical extended mode profile is shown in Fig. [4j plotted 
both as a function of the physical coordinate x for large 
supercell and also in the unit cell of the superspace (in- 
set). In the inset superspace plot, one can clearly see 
the predicted field oscillations perpendicular to the slice 
plane, as well as a slower oscillation rate (inversely pro- 
portional to the frequency) parallel to the slice. In the 
plot versus x, one can see the longer-range quasi-periodic 
structure that arises from how the slice wraps around the 
unit cell in the superspace. The factor of three to four 
long-range variations in the field amplitude are sugges- 
tive of the critically localized states (power -law d ecay) 
that one expects to see in such quasicrystals.'^'^^l^ 



B. Defect modes 

Much of the interest in quasicrystal band gaps, similar 
to the analogous case of band gaps in periodic structures, 
centers around the possibility of localized states: by in- 




FIG. 5: Dielectric for the Fibonacci chain with e — 2.56 (light 
blue), and a defect — an additional e = 8.0 layer, shown in 
gray. 



troducing a defect in the structure, e.g. by changing the 
thickness of a single layer, one can create exponentially 
localized states in the gap.^^^ In periodic systems, be- 
cause such defects break the periodicity, they necessitate 
a larger computational cell, or supercell, that contains 
many unit cells. In quasicrystal systems, however, one 
can introduce a localized defect without changing the 
higher-dimensional periodicity, and therefore compute lo- 
calized defect modes with the same superspace method 
and computational cell. 

Ideally, if one had infinite spatial resolution, a defect 
in the crystal would be introduced as a very thin pertur- 
bation parallel to the slice direction. As the thickness of 
this perturbation goes to zero, it intersects the physical 
slice at greater and greater intervals in the physical space, 
corresponding to localized defects that are separated by 
arbitrarily large distances. In practice, of course, the 
thickness of the perturbation is limited by the spatial 
resolution, but one can still obtain defects that are very 
widely separated — since the associated defect modes are 
exponentially localized, the coupling between the defects 
is negligible. In other words, one effectively has a very 
large supercell calculation, but expressed in only the unit 
cell of the higher-dimensional lattice. 

As an example, we changed an £ = 2.56 layer to 
e = Ed one place in the Fibonacci quasicrystal. The 
corresponding superspace dielectric function is shown in 
Fig. [sj where the defect is introduced as a thin (0.02a) 
strip of Ed parallel to the slice direction. We compute 
the band structure as a function of the defect dielectric 
constant As = Sd — 2.56, varying it from the normal di- 
electric Ed = 2.56 up to Ed = 11. The thickness of the 
defect in the unphysical direction was fixed to be 0.02. 
The resulting eigenvalues as a function of Ae are shown 
in Fig. |6] for two different spatial resolutions of 50 (blue) 
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FIG. 6: Varying the defect epsilon for resolutions 50 (blue) 
and 100 (red). The thickness of the defect is fixed to 0.02 
lattice constants. The number of spurious modes increases 
with the resolution, the true defect state being the lowest of 
these modes. 

and 100 (red) pixels/a. When the resolution is 50 the 
defect is only one pixel thick, the discretization effects 
might be expected to be large, although the frequency is 
within about 2% of the higher-resolution calculation. At 
the higher resolution, the frequency of the mode is con- 
verging (it is within 0.3% of a resolution-200 calculation, 
not shown). However, at the higher resolution there is 
a second, spurious mode due to the finite thickness (2 
pixels) of the defect layer — this spurious mode is easily 
identified when the field is plotted Fig. [Tj^bottom), be- 
cause it has a sign oscillation perpendicular to the slice 
(which would be disallowed if we could make the slice 
infinitesimally thin). 

The defect modes for the resolution 100 are plotted in 
Fig. [7| for both the real and the spurious modes, versus 
the physical coordinate (x) and also in the superspace 
unit cell (insets). When plotted versus the physical co- 
ordinate X on a semilog scale, we see that the modes are 
exponentially localized as expected. The defect mode ap- 
pears at multiple x values (every ~ 20a on average) be- 
cause the defect has a finite thickness — the physical slice 
intersects it infinitely many times (quasiperiodically), as 
discussed above. The spurious mode (bottom panel) is 
also exponentially localized; it has a sign oscillation per- 
pendicular to the slice direction (inset) which causes it 
to have additional phase differences between the different 
defects. 

The advantages of the higher-dimensional (superspace) 
calculation over a traditional supercell calculation are 
more tenuous for this sort of defect calculation, because 
the exponential localization means that a relatively small 
supercell can be employed. On the other hand, this is an 
illustration of the versatility of the superspace approach 
and is a powerful tool for studying quasiperiodic struc- 




FIG. 7: Semi-log plots of the magnetic field magnitude Hz for 
the lowest (top) and highest (bottom) defect state for the con- 
figuration shown in Fig. [5] Insets: Two-dimensional super- 
space visualizations of the defect states. Note the additional 
node in the lower figure (corresponding to an unphysical os- 
cillation) . 

tures with or without defects. 



C. Continuously varying the cut angle 

The cut-and-project construction of quasicrystals pro- 
vides a natural way to parameterize a family of periodic 
and quasiperiodic structures, via the cut angle (j). It is 
interesting to observe how the spectrum and gaps then 
vary with (f). 

As (j) is varied continuously from 0° to 45°, the struc- 
tures vary from period a to quasi-periodic lattices (for 
tan (j) irrational) to long-period structures (tan cj) rational 
with a large denominator) to a period a\/2 crystal. As 
we change (/>, we rotate the objects in the unit cell, so 
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FIG. 8: Projected band structure vs. cut angle 0, showing 
different one-dimensional quasicrystal realizations. 



that they are always extruded along the y direction with 
a length equal to the projection of the unit cell onto y 
[a(sin 0+cos ^)], corresponding the usual cut-and-project 
construction.^^ In this case, the spectrum varies continu- 
ously with ^, where the rational tan (j) correspond to "ra- 
tional approximants" of the nearby irrational tan 0.1^3211 
For a general unit cell with a rational tan 0, the physical 
spectrum might depend on the slice offset y and hence 
different from the total superspace spectrum, but this is 
not the case for dielectric structures like the one here, 
which satisfy a "closeness" condition^'^ (the edges of the 
dielectric rods overlap when projected onto the Y direc- 
tion). This makes the structure y-independent even for 
rational slices^^ The resulting structures are shown in the 
bottom panel of Fig. [S] for three values of (j). 

The corresponding photonic band gaps are shown in 
the top panel of Fig. [Sj as a continuous function of (j). 
Only the largest gaps are shown, of course, since we are 
unable to resolve the fractal structure to arbitrary res- 
olution. As might be expected, there are isolated large 
gaps at (/) = 0° and (j) = 45° corresponding to the simple 
ABAB • • • periodic structures at those angles (with pe- 
riod a and a/V^, respectively, the latter resulting from 
two layers per unit cell). The (j) = 45° gap is at a higher 
frequency because of its shorter period, but interestingly 
it is not continuously connected to the = 0° gap. The 
reason for this seems to be that the two gaps are dom- 
inated by different superspace reciprocal lattice vectors: 



(1,0) • 27r/a for = 0°, and (1,1) • 27r/a for = 45°. 
For intermediate angles, a number of smaller gaps open 
and then close. If we were able to show the spectrum 
with higher resolution, we would expect to see increas- 
ing numbers of these smaller gaps opening, leading to 
the well-known fractal structure that arises e.g. for the 
Fibonacci crystal. 

This variation of gaps as a function of (f) provides in- 
teresting possibilities for band-gap engineering. First, 
we see that we can get large gaps that are close to 
one another in the spectrum, whereas in a typical one- 
dimensional quarter-wave stack the gaps are at integer 
multiples of a given frequency. Even more complex 
combinations of gaps may arise for higher dimensional 
superspaces, since by including additional incommensu- 
rate reciprocal lattice vectors one may generate addi- 
tional nearby gaps. Yet another interesting possibility 
would involve optimizing the layer thicknesses as a func- 
tion of (j) so as to maximize the largest band gap at every 
(f) (or some other objective). The dielectric layers cho- 
sen in Fig. [s] were of fixed thickness 0.37a (and variable 
height equal to the projection of the unit cell along the 
y direction as described above). 



V. CONCLUDING REMARKS 

We have presented a numerical approach to computing 
the spectra of photonic quasicrystals by directly solving 
Maxwell's equations extended to a periodic unit cell in 
higher dimensions, allowing us to exploit Bloch's theorem 
and other attractive properties of computations for peri- 
odic structures. In doing so, we extended the conceptual 
approach of cut-and-project techniques, which were de- 
veloped as a way to construct quasicrystals, into a way 
to simulate quasicrystals. Compared to traditional su- 
percell techniques, this allows us to capture the entire 
infinite aperiodic quasicrystal in a single finite compu- 
tational cell, albeit at only a finite resolution. In this 
way, the single convergence parameter of spatial resolu- 
tion replaces the combination of resolution and supercell 
size in traditional calculations, in some sense uniformly 
sampling the infinite quasicrystal. The resulting com- 
putations, applied to the test case of a Fibonacci qua- 
sicrystal, display the unique features of quasicrystals in 
an unusual fashion, in terms of higher-dimensional band 
structures and visualization techniques. This technique 
also allows defects and variation of cut angle (continu- 
ously varying between periodic and aperiodic structures) 
in a straightforward way. 

In future work, we plan to apply this approach to mod- 
eling higher-dimensional quasicrystal structures, where 
computing the spectrum is currently more challenging 
using existing supercell techniques. To make a higher- 
dimensional superspace calculati on pra ctical, one must 
use iterative eigensolver methods^^^ rather than the 
simple dense-matrix techniques employed for our test 
case. Iterative techniques are most efficient for comput- 
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ing a few eigenvalues at a time, and so it will be useful 
to employ iterati ve ni ethods designed to compute "inte- 
rior" eigenvalues, S^l^ allowing one to search directly for 
large gaps without computing the lower- lying modes. Al- 
ternatively, numerical techniques have been developed, 
based on filter-diagonalization methods, to directly ex- 
tract the spectrum of many eigenvalues without comput- 
ing the corresponding eigenvectors.^^ 

Appendix 

In this appendix, we give an explicit derivation of the 
fact that an "irrational" slice densely fills the superspace 
unit cell, or rather a definition of the necessary conditions 
to be an "irrational" slice. These concepts are widely 
used in the quasicrystal literature, but a precise definition 
seems hard to find (one commonly requires that all of 
the Miller indices have incommensurate ratios, but this 
condition is stronger than necessary). 

Without loss of generality, we can consider the unit cell 
in the superspace Z = R"^ to be the unit cube (related to 
any lattice by an affine transformation) with lattice vec- 
tors along the coordinate directions. The physical slice X 
is (i-dimensional, and it will be convenient to write the co- 
ordinates of a vector z as z = (si, . . . , 5^^, ti, . . . , tn-d) = 
(s, t). By taking every coordinate modulo 1, we can map 
X to a set X consisting of X's intersection with each unit 
cell. We wish to show necessary and sufficient conditions 
for X to densely fill the unit cell. 




FIG. 9: An example two-dimensional superspace slice X (left) 
and the projected slice modulo 1 into the unit cell X in the 
s (middle) and t directions (right), along with the intersec- 
tions T (middle, red) and T (right, red) of X with the s = 
hyperplane. 



Assuming that the slice is not orthogonal to any of 
the coordinate axes (as otherwise it would clearly not 
densely fill the unit cell), we can parameterize the points 
z of X so that the last n — d coordinates (ti, . . . , tn-d) 
are written as a linear function t(si, . . .Sd) = t(s) of the 
first d coordinates. 



Consider the set T in W^~^ formed by the t(s) coor- 
dinates of X when the components of s take on integer 
values. This is a subset of X, and the corresponding set 
T formed by taking t G T modulo 1 is a subset of X. The 
key fact is that X is dense in the n-dimensional unit cell 
if and only if T is dense in the (n — (i) -dimensional unit 
cell, and this is the case that we will analyze. This equiv- 
alence follows from the fact that X is simply T translated 
continuously along the slice directions (every point in X 
is related to a point in T by a simple projection). The 
set T is a lattice in W^~^ consisting of all integer linear 
combinations of the basis vectors tk = t{sj = Sj^), since 
t(s) is a linear function. 

For each basis vector t/c, it is a well-known factP^that 
if it consists of m incommensurate irrational components, 
the set of integer multiples it^ modulo 1 will densely fill 
an m-dimensional slice of the unit cell. More precisely, 
write tk = Z)j,=i...mfc ^k^l + ^fc' where the and qk 
have purely rational components and the {aj} are in- 
commensurate irrational numbers, and m/e is therefore 
the number of incommensurate irrational components of 
t/e. Then the set of integer multiples of modulo 1 
densely fills an m/g-dimensional slice of the unit cell of 
M^~^. The basis vectors of this slice are precisely the 
vectors b-^, which are rational and therefore commensu- 
rate with the basis vectors of R^~^, while the vector q/. 
is simply a rational shift. This slice therefore cuts the 
unit cell of R^~^ a finite number of times. 



The set T is then obtained as the direct sum of these 
dense slices for all n — d vectors t/^. This is then dense 
if and only if the set of vectors {b"^}i=i"""^'' spans R^"'^. 
In other words, an "irrational slice," which densely fills 
the unit cell, is one in which there are n — d independent 
incommensurate slice components as defined above. 
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